On the behavior of a 2D Heisenberg antiferromagnet at very low temperatures 
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We present an analytical result for the ratio of the physical correlation length in a 2D Heisenberg 
antiferromagnet on a square lattice, and the one which is actually computed in numerical simulations. 
This last correlation length is deduced from the second moment of the structure factor at the 
antiferromagnetic momentum Q. We show that the ratio is very close to one in agreement with 
previously obtained numerical result of the 1/iV expansion. 



The two-dimensional Heisenberg antiferromagnet on 
a square lattice is one of the most extensively studied 
systems in condensed-matter physics. The interest to 
this model is two-fold. On one hand, Heisenberg an- 
tiferromagnet models a large number of real materials 
including parent compounds of high-T c superconductors. 
On the other hand, its low-energy physics is adequately 
described by a field-theoretical a— model thus allowing 
one to find similarities between condensed-matter physics 
and field theory. 

The low-temperature behavior of the Heisenberg an- 
tiferromagnet is understood in great detail. For short- 
range interaction, the ground state is ordered unless one 
fine tunes the couplings between nearest and further 
neighbors. The ordered ground state is characterized by 
a sublattice order parameter, N$, spin stiffness, p s , and 
transverse susceptibility, Xl = c ~ 2 Ps, where c is the spin- 
wave velocity. At any finite temperature, however, the 
system is disordered due to thermal fluctuations. The 
disordering means that equal-time spin-spin correlation 
function decays exponentially with the distance, as e~ r ^ . 
The length scale £ is the physical spin correlation length. 
Various approaches to 2D antiferromagnets all predict 
that in the renormalized-classical region (T <§; p s ) which 
we only consider here, £ is exponentially large in T at 
low T and behaves as £ ~ ex P(^r^)- Equal-time 
spin correlations at large distances can also be described 
by a static structure factor S(k) for k near the antiferro- 
magnetic momentum Q = (jr, it). In the disordered spin 
state, S(Q) scales as £ 2 and is therefore also exponential 
in T. 

The exponential temperature dependences of £ and of 
S(Q) have been verified in numerical simulations and 
by analyzing the neutron scattering and NMR data for 
La,2CuOi and .S^CwC^C^- The accuracy of numer- 
ical simulations is so high, however, that one can not 
only check the temperature dependences but also com- 
pare the absolute values of the structure factor and the 
spin correlation length. There is, however, one subtlety 
in the numerical analyses - in numerical simulations, one 
measures not the physical spin correlation length, but 
another length scale which is also exponentially small 



in T, but generally differs from £ by a constant factor. 
Indeed, the physical spin correlation length can be ex- 
tracted from the form of S(k) simply because Fourier- 
transform of S(k) yields the spin correlation function in 
real space. However, one can easily make sure that the 
long-distance behavior of S(r) is associated with the form 
of S(k) along the imaginary fc-axis. Moreover, the inverse 
spin correlation length is the scale at which S(k) has a 
pole for imaginary k: S~ 1 (k = iS,^ 1 ) = 0. In numerical 
simulations, however, the structure factor is evaluated 
only for real values of momentum k. Accordingly, a dif- 
ferent definition of the correlation length is employed - it 
is identified as a second moment of S(k) for k = Q, i.e., 
as£ = (-S- 1 (Q)dS(k)/dk 2 \ k ^ Q ) 1 / 2 . For the Lorentzian 
form of S(k), S(k) cx ((Q - kf + m 2 ,)" 1 , both £ and £ 
are equal to the mean-field spin excitation gap mo and 
are therefore undistinguishable. However, the 1/N calcu- 
lations for the O(N) a— model rigorously demonstrated 
that S(k) has a Lorentzian form only in the limit N — > oo 
while for arbitrary N, and, in particular, for physical 
A = 3, S(k) possesses a much more complex dependence 
on k which arizes due to temperature-independent 1/N 
corrections. In this situation, £ and the physical spin 
correlation length £ differ by some constant factor. This 
factor appears to be a relevant one as several groups 
recently performed a detailed comparison of £, calculated 
numerically at very low T, with the exact expression for £ 
obtained some time ago by Hasenfratz and Niedermayer 
(see below). In this analysis, they assumed that £ and 
£ are almost the same. A priori, however, there are no 
reasons for such coincidence. 

We start by quoting the seminal result M that 
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This result is based on numerical results for N — 3 and 
N = 4 @ and on 1/N expansion for 0(A) cr-model [[§§. 
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For N = 3, this yields ^/m = (8/e) w 2.94. No such 
exact expression, however, exists for £. The 1/iV ex- 
pansion for O(N) a model (which in the renormalized- 
classical region holds in powers of 1/(N — 2) yields 



r 1 = r 1 (i + o.oo3/(iv-2)), 



(3) 



where factor 0.003 arises from numerical evaluation of 
some complicated integrals. A formal application of this 
result to a physical case of N = 3 yields almost identical 
values for £ and £. It is known, however, that a special 
care has to be taken in this procedure as not all terms 
which appear to first order in 1/(N — 2) actually con- 
tribute at N = 3. To illustrate this point, consider large 
iV expansion for An analytical evaluation of the 

first 1/N correction yields 
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where is the Euler constant. Comparing this formula 
with the exact result, Eq. ([!]), we notice that the Euler 
constant accounts for the appearance of the T function in 
(§: r(l + l/(iV-2)) = l- lE /{N-2) + 0{l/(N-2) 2 ). 
As F(2) = T(l) = 1, the term with the Euler constant 
does not contribute to £ _1 /m for the physical case of 
N = 'i. The danger is that the same might also happen 
for the rescaling factor between £ and £, i.e., that the 
ratio obtained numerically to first order in 1/N may in 
fact contain the Euler constant which would mask the 
actual value of the ratio. 

In the present communication we address this issue. 
We compute explicitly the T— independent 1/N correc- 
tions to both correlation lengths and show that the rescal- 
ing factor does not contain the Euler constant. This im- 
plies that the 1/N result for the ratio is very likely to be 
trustable, and the rescaling factor between two correla- 
tion length is very close to one. 

As we already said, in our consideration, we will 
use heavily results of the 1/N expansion for O(N) 
a— model 0,^]. The point of departure for this analy- 
sis is the mean- field theory which is exact at N — oo. At 
the mean-field level, the static structure factor is given 
by 
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and m = (T/c)e- 2 *P°/( N7 "> . Consider now finite N. It 
has been shown in |^,|| that there exist two different 
types of 1/(N — 2) corrections: singular ones which con- 
tain log(T/mo) and log(log(T/mo)), and regular ones, 
which at small T account for the T-independent renor- 
malizations of S(k) and m. There is a large amount of 
confidence that logarithmical and double logarithmical 
perturbation series are geometrical and therefore can be 
simply exponentiated M. Regular l/(N—2) corrections 



require a special care, as we just demonstrated. Col- 
lecting both singular and regular 1/iV corrections and 
exponentiating the singular ones, one obtains [H 



S(k) = 2tt7V 2 
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N-2 V 2tt Ps 
where for k comparable to the inverse correlation length 
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P(k) 



Zm 2 + (Q-fc) 2 + S(fc)' 
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Here m is given by (0) and Z and S(fc) oc (Q—k) 2 account 
for a temperature independent 1/(N— 2) corrections. For 
N = 3, the functional forms of S(k) and m obtained this 
way fully agree with the ones obtained in the perturbative 
RG approach ||. 

The expressions for Z and S(fe) have been obtained 
in pj but not explicitly presented in that paper. Here we 
list the catalog of the results which we will need: 
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For the two correlation lengths we then obtain 

S(fc = im) . 



(8) 



C = m 2 Z(l 
C 2 = m 2 Z(l 



E(fc) 
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Performing simple manipulations, we find 
r 1 =m(l+^(21o g 2-l-A) 

where 



(9) 



(10) 



, , . . x + v / ?+4 ( (x- Vx 2 + 4) 2 
A= I dxloglog 



2y/x 2 + 4 



(11) 



Introducing t — log x +^^ 2 + 4 and integrating by parts, 
we immediately obtain A = — (j E + log 2). Substituting 
this result into (RToT) we recover Eq. (H). For the ratio of 
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£ and £, the same manipulations yield a more complex 
expression 



where 7 = 37 2 - 14J 3 + I8J4, and 



(12) 



I„ 



dt- 



sinh t log £ 
(2 cosh i- l) n 



Notice that all integrals /„ are convergent. Physically, 
this implies that T— independent factors in both correla- 
tion lengths are determined by system behavior at energy 
scales which are much smaller than the upper cutoff. At 
this scales, the behavior is universal and therefore overall 
factors in £ and £ are also universal numbers. Note by 
passing that similar integrals appear in the calculations 
of the d.c. Hall conductivity near fractional quantum 
Hall critical point |J. 

To evaluate the integrals (|l3|), we introduce auxiliary 
function 3> n (t) = sinhi log 2 (— t)/(2cosht — 1)™ and inte- 
grate over a contour which consists of a circle of an 
infinite radius and a cut along positive real t. The con- 
tour integral yields — 47ri/„ and simultaneously it is equal 
to the sum of the residues (modulo 2ni) of the poles along 
imaginary t— axis. Performing calculations and making 
use of the summation formula 



log(27rn + 7r(l - a)) log(27ro + 7r(l + a)) 
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we can explicitly pull out the Euler constant from the 
integrals 
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where 



h = 
h = 
h = 
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The analytical expressions for the integrals 1%^ have 
been obtained very recently JTo| |, whereas I2 was known 
for some time U]. It turns out that these integrals 
can be expressed in terms of the derivatives of the 
Hurvitz Zeta function d(,{x,a) / dx at x = 0,-1,-2 and 
a = 1/6,1/3,1/2,2/3,5/6,1. For x = 0, d((x,a)/dx = 



log(r(a)/v27r). By analogy we can introduce gen- 
eralized Gamma functions via e?£(x, a)/dx\ x= - n = 
log(r_ n (a)/v / 27r)- Using this definition, we obtain a 
symbolic representation 

h = -logi? ; ^3 = 61ogi?_i - 3 log 6; 

7 4 = 361ogi?_ 2 +31og6 (17) 



(13) where 
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The expression for R and can be reduced to 

an explicit closed-form expression containing, e.g., log- 
arithms of Gamma- functions, but no Eulcr constant |l2| . 
No closed-form expression in terms of tabulated func- 
tions is known for i?_2- It is however quite reasonable to 
assume that it also does not contain Euler constant. 

Assembling now all contributions to I and substituting 
the result into ([[2]) we find that the Euler constant is 
cancelled out. The rest is combined into 
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This expression is the central result of the paper. 

The next issue is how to account for the higher-order 
terms in 1/(N — 2). Here we use the same assumption as 
was proven to work for £, namely that after the Euler con- 
stant is subtracted, the rest of the regular 0(1/ (N — 2)) 
correction can be exponentiated. Using this assumption, 
we finally obtain 



f 2 = e {&Rl /3 RZlRzT) 1/iN - 2) 



(20) 



For TV = 3 this yields f 2 /£ 2 = 0.993, i.e. the ratio is 
indeed very close to one. This result may sound intu- 
itively obvious, but is not based on any apparent phys- 
ical reasons and therefore had to be verified by explicit 
calculations. This is what we did. 

The extreme closeness of the ratio £/£ to 1 is consis- 
tent with recent claims that at T — > 0, the numerically 
computed spin correlation length J^,|| approaches the 
Hasenfratz-Niedermayer result, Eq.(|l]). 

Using our expressions for Z and £, we can also com- 
pute the overall factor for the structure factor S(Q) which 
is also studied in numerical simulations. The numer- 
ical evaluation of the first 1/(N — 2) correction yields 
P(Q) = l/(Zm 2 ) = £ 2 (1 + 0.188/(JV - 2)). Exponenti- 
ating this result, one obtains P(Q)/^ 2 ~ 1.2. Numerical 
simulations Q , on the other hand, reported that the ac- 
tual rescaling factor is more than three times larger than 
this number. We performed analytical calculations along 
the same lines as above and obtained 
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where I2 is given by (|l5|). Substituting the result for I 2 
into this expression, we obtain 



P(Q) = C 2 [i 
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(- 7B + 31og(i? )-log2)] (22) 



This result coincides with the one obtained earlier by 
Campostrini and Rossi ||, and cited previously in [ jl3| . 
We see now that Euler constant is present in the pertur- 
bation series, i.e. one cannot simply exponentiate the 
lowest-order result. Using the same procedure as be- 
fore, i.e., treating 7^ as coming from the expansion of 
r(l + 1/(N — 2)), and exponentiating the rest of (|22]), we 
obtain 



p(Q) = 2 1 /( 2 - Ar ) r( 
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For N = 3, this yields P{Q)^ 2 = 2.149 which is about 
two times larger than the result obtained by formally ex- 
ponentiating the whole 1/(N — 2) correction. Still, how- 
ever, this result does not fully agree with quantum Monte 
Carlo simulations at low T which reported P (Q)£,~ 2 ~ 4 
for both S = 1/2 [§ and S = 1 @. The series expansion 
results [|U reported somewhat smaller P(Q)£~ 2 ~ 3.2 
for S = 1/2. The reason for the discrepancy is not clear 
to us. Possibly numerical simulations for S(Q) were per- 
formed not deep inside the asymptotic scaling regime at 
T — > 0. Another possibility is that something may be 
wrong with the exponentiation of the first l/N correction 
to S(Q), though this is unlikely in view of the results for 
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